from datetime import datetime

import pandas as pd
import numpy as np
import xarray as xr

import ttide
from ttide.t_tide import t_tide
from scipy.optimize import curve_fit


from ppgnss import gnss_time
from ppgnss import gnss_utils
from ppgnss import gnss_geodesy

from plot_baseline import read_seal_file, read_gamit_solution
import matplotlib.pyplot as plt
font_path = '/Library/Fonts/PingFang.ttc'  # 请替换为实际的字体文件路径

plt.rcParams['font.sans-serif'] = ['Arial Unicode MS']
# plt.rcParams['font.family'] = 'Songti SC'
# plt.rcParams['font.sans-serif'] = ['Songti SC']
# plt.rcParams['font.sans-serif']=['SimHei','Songti SC','STFangsong']

gamit_solution = "../baseline/solution.dat"

zihe_excel_file_path_02 = "../data/cczihe02.xlsx"
zihe_excel_file_path_03 = "../data/cczihe03.xlsx"
excel_file_path = zihe_excel_file_path_02

seal_down_file = "../data/sea_levels/downstream_7308.txt"
seal_up_file = "../data/sea_levels/upstream_3020.txt"

lat1, lon1 = 30.22863399, 120.73933035
lat2, lon2 = 30.22401975, 120.74423075

x1, y1, z1 = gnss_geodesy.blh2xyz(lat1, lon1, 0)
x2, y2, z2 = gnss_geodesy.blh2xyz(lat2, lon2, 0)

north, earth, up = gnss_geodesy.dxyz2neu([x2-x1, y2-y1, z2-z1], [x1, y1, z1])
alpha = np.arctan(north/earth)

df = pd.read_excel(excel_file_path)
# xr_gg = read_gamit_solution(gamit_solution, 2022, 118)
xr_seal_down = read_seal_file(seal_down_file)

xr_seal_up = read_seal_file(seal_up_file)
xr_seal_down, xr_seal_up = xr.align(xr_seal_down, xr_seal_up, join="inner")
time_from, time_to = "2023-01-29 00:00:00", "2023-03-09 00:00:00" # 50 days
xr_delta = xr_seal_up.loc[time_from:time_to] - xr_seal_down.loc[time_from:time_to]
width, height = gnss_utils.cm2inch(18), gnss_utils.cm2inch(10)
fig, axes = plt.subplots(nrows=3, ncols=1, figsize=(width, height), sharex=True, sharey=True)

plt.subplots_adjust(left=0.08, right=0.98, bottom=0.12, top=0.98, hspace=0)

coord_time = xr_seal_down.loc[time_from:time_to].coords["time"]
xind = [(_-coord_time[0])/np.timedelta64(1, "D") for _ in coord_time]
linewidth=1
h0, =axes[0].plot(xind, xr_seal_down.loc[time_from:time_to], '-', linewidth=linewidth, label=u"闸外水位")
axes[0].set_ylabel(u"闸外水位 (m)")
h1, =axes[1].plot(xind, xr_seal_up.loc[time_from:time_to], '-',linewidth=linewidth, label=u"闸内水位")
axes[1].set_ylabel(u"闸内水位 (m)")
# xind = xr_delta.loc[time_from:time_to].coords["time"]
h2, =axes[2].plot(xind, xr_delta.loc[time_from:time_to], "-",linewidth=linewidth, label=u"内外位差")
axes[2].set_ylabel(u"内外位差 (m)")
axes[0].grid(axis='y', linestyle='-')
axes[1].grid(axis='y', linestyle='-')
axes[2].grid(axis='y', linestyle='-')

# xind = xr_seal_up.loc[time_from:time_to].coords["time"]
# handles, labels = plt.gca().get_legend_handles_labels()

xticks = range(0, 50, 5)
xlbls = []
for idy in xticks:
    start_time = coord_time[0].values
    current_time = start_time + np.timedelta64(idy, "D")
    # print(current_time)
    lbl = str(current_time)[5:10]
    xlbls.append(lbl)
yticks = np.arange(-2.5, 7.6, 2.5)
ytklbls = [f"{_:.1f}" for _ in yticks]
plt.xticks(xticks, xlbls, fontsize=10)
plt.yticks(yticks, ytklbls[:-1]+[None], fontsize=10)
plt.xlim(xind[0], xind[-1])
plt.ylim((-2.5, 7.5))
# plt.legend([h0, h1, h2], [u"闸外水位", u"闸内水位", u"内外位差"], loc='upper center', ncol=3, fontsize=10,bbox_to_anchor=(0.45, 1.26))
plt.xlabel(u"时间 (2023年, 月-日)", fontsize=10)

fig_filename = "figures/sea_level.png"
plt.savefig(fig_filename, dpi=300)
plt.show()

data = {"down": xr_seal_down,
        "up": xr_seal_up,
        "delta": xr_delta}

filename = "sea_level_0129_0309.npo"
gnss_utils.saveobject(data, filename)
print(filename + " is saved!")